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In full flux-surface computer studies of tokamak edge tur- 
bulence, a spectrum of shear flows is found to control the 
turbulence level and not just the conventional (0,0)-mode 
flows. Flux tube domains too small for the large poloidal 
scale lengths of the continuous spectrum tend to overesti- 
mate the flows, and thus underestimate the transport. It is 
shown analytically and numerically that under certain condi- 
tions dominant (0,0)-mode flows independent of the domain 
size develop, essentially through Bose-Einstein condensation 
of the shear flows. 
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The energy confinement of tokamaks is mainly con- 
trolled by small-scale (~cm) turbulence giving rise to 
the "anomalous transport" . Analytical and computer- 
aided studies have found that the anomalous transport 
in turn is often controlled by fluctuating "zonal flows" 
jl],D, poloidal shear flows, which are assumed to have 
zero poloidal and toroidal mode number, but have ra- 
dial scales similar to the turbulence. The present paper 
deals with the question what happens to the shear flows 
when the turbulence scale-lengths become very small 
compared to the plasma size such as in the tokamak edge, 
in particular in future large machines. In this limit, the 
shear flows can either contain a finite (0, 0)-mode com- 
ponent, or may loose their global character and change 
into vortices with finite poloidal scale length, as will be 
demonstrated by numerical full flux-surface edge turbu- 
lence studies. Regarding their large poloidal and parallel 
but small (similar to the turbulence) radial scale-lengths, 
these vortices should not be regarded as drift-waves or 
convective cells || but rather as poloidally localized shear 
flows. 

For cost reasons, the domains of turbulence simula- 
tions are usually thin flux tubes [Q or tokamak sectors 
f|, equivalent to flux tubes with special boundary con- 
ditions. The flux tube dimensions perpendicular to the 
magnetic field are of the order of ~ 10cm and they extend 
several ~ 10m along the magnetic field to accomodate 
the prevalent turbulent structures. For poloidally local- 
ized shear flows, however, these computational domains 
are not adequate and the flows always appear to extend 
across the complete flux tube. Since the (0, 0) mode is 
not damped as the other modes, it may therefore exert 
a strong stabilizing effect on the turbulence, which un- 
derestimates the transport compared to a true full flux- 
surface simulation. For the core turbulence, which has 
relatively large scales, full torus simulations |^[7j] exhibit 
zonal flows extending over the complete flux surfaces. 
However, even for these scenarios it is not clear, whether 



this remains true for much larger ratio of flux surface 
circumference to turbulence scale length or whether the 
flows have a finite scale length in poloidal direction. 

Following the numerical results, an analytic model for 
the shear flows is described, in which their poloidal and 
radial wavenumber spectra are controlled by the interplay 
of damping by the collisional electron response and ion 
dissipation, the linear response of the turbulence to the 
flows, and the excitation of flows by random fluctuations. 
Under certain conditions, a non-zero fraction of the flow 
energy is generated as (0, 0)-mode flows, regardless of the 
system size. The mechanism is analogous to the Bose- 
Einstein condensation (BEC). The three effects acting 
on the flows take the role of absorption, stimulated, and 
spontaneous emission. For the BEC, a macroscopic frac- 
tion of the quanta is eventually scattered into the ground 
state because the state density near the ground state is 
too low to hold sufficiently many quanta under the preva- 
lent conditions. For the flow system, the turbulence and 
shear flows form a feed back loop, which regulates the 
shear flow energy to the level needed for turbulence sat- 
uration. Condensation into the (0, 0) flow component 
occurs when a threshold in required flow energy is ex- 
ceeded, and the ra^O modes are unable to hold it. 

Numerical results. - We discuss the results of tur- 
bulence simulations of the three dimensional electrostatic 
drift Braginskii equations with isothermal electrons (a 
subset of the equations of Ref. || ) for two different cases: 

(a) the predominant instability is the resistive balloon- 
ing mode with the nondimcnsional parameters ay = 0.2, 
e n = 0.08, q = 5, r = 1, rji = 1, § = 1; (b) there is a 
significant contribution from ITG modes with ad = 0.4, 
rji = 3 and the other parameters as in (a). The radial 
domain width in terms of the resistive ballooning scale 
length, Lo, in (a) was 2ALq and in (b) 48Loj the width 
Lg perpendicular to r and B was 24Lo [only for (a)], 
192Lo, 384Lo, and 768Lo (the corresponding tokamak 
minor radius is a = Lgq/(2ir)). For a definition of these 
parameters and units see Refs. The parameters of 
the largest domain are consistent with the physical pa- 
rameters R — 3 m, a — 1.5 m, L n — 12 cm, qo = 3.2, 
n = 3.5 x 10 19 m" 3 , Z cfi = 4, B = 3.5 T, and for (a) 
T = 100 eV, L = 5.1 mm, p s = 0.58 mm and for (b) 
T = 200 eV, L = 3.6 mm, p s = 0.82 mm. The perpen- 
dicular grid step size was A = 0.19Lo (a) and A — 0.38Lo 

(b) . Parallel to the magnetic field 12 points per poloidal 
connection length were sufficient due to the large parallel 
scales of the ballooning modes. The largest runs had a 
grid of 128 x 4096 x 12. 

The dependence of the average (0, 0) shear flow energy 
density on the domain size, Lg, is compared for the two 
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cases in Fig. [l| In contrast to case (b), the shear flows 
in (a) are apparently not condensed into the (0, 0) mode 
since its energy density decreases proportional to 1/Lg oc 
1 /a, as is expected when a given shear flow energy density 
is distributed equally among an increasingly dense set of 
modes. 

The kg spectrum of the flow velocity, v = vg — d r <j>, 
for the Lg = 768L runs [for case (a) see Fig. || exhibits 
a rise at low kg associated with the shear flows, different 
from the microturbulence fluctuations at kg ~ 1. The 
square amplitude of the m = mode in case (a) and (b) 
is 0.3 and 7 times, respectively, the total shear flow am- 
plitude, suggesting strong condensation for (b). In both 
cases, the typical poloidal scale length of the to ^ shear 
flows is roughly a factor 10 greater than the scales of the 
turbulence. Failure of the computational domain to ac- 
comodate the scales of the uncondensed shear flows in 
case (a) results in an overestimate of the shear flow am- 
plitude, and hence in an underestimate of the anomalous 
transport. The particle flux for Lg = L r = 2ALq was 
found to be 25% lower than for Lg = 768L . 

Analytic model. - As the first ingredient of a quali- 
tative model for the poloidal shear flow spectra, we cal- 
culate the linear dispersion relation for finitely elongated 
shear flows. For clarity, in the linear electrostatic vortic- 
ity equation (with the plasma parameters absorbed into 
the units, see, e.g., M), 



vi(a t + 7 )^ + af0 = o, 



(i) 



we neglect temperature fluctuations, parallel ion velocity, 
drift effects, curvature and magnetic fluctuations. These 
effects can lead to a real frequency (e.g., geodesic acous- 
tic modes [|lO|Jllf| ) and to a coupling to parallel sound 
waves or Alfven waves. The dissipative effects are the 
flow damping 7 due to the ion dissipation assumed inde- 
pendent of the wavenumber and the damping due to the 
resistive electron response. As we will see below, for a po- 
tential condensation of the shear flows into global modes 
only a small region around a certain radial wavenumber, 
ko, is important, which we set to one in (|l|) since its abso- 
lute value is not important. Because of the large poloidal 
wavelengths of the flows, we approximate — V\ w k$ = 1 
and obtain the dispersion relation 



i ( 7 + ki 



q(r) 



(2) 



The damping by the parallel resistive electron response 
is weak if either m = n = holds, or r is near a resonant 
surface defined by m — nq(r mn ) = 0. Focusing on a thin 
region around r — ro we obtain k\\ ~ mao(r — r mn ), 
cxq = —q'(ro)/q{ro) 2 - Hence the resistive flow damping is 
proportional to to 2 , which is the reason for the poloidal 
elongation of the flows, i.e., their low mode numbers. 

As reaction to a shear flow p|,p| Jl2|jr^ ] the microturbu- 
lence may in turn influence the flows via the Reynolds 
stress [|l4 15 or the Stringer- Winsor mechanism due to 



poloidal pressure asymmetries |I^,|l6[. Restricting our- 
selves to linear response theory, we assume a (coherent) 
flow amplification rate g(k r ) depending only on the radial 
wavenumber k r , because of the large poloidal correlation 
lengths of the shear flows. With the (incoherent) random 
forcing, /, representing the effect of the turbulence fluc- 
tuations, the equation for the flow amplitude in frequency 
space has the form of a Langevin equation, 



d t v = -icuhnv + g{k r )v + f. 



(3) 



From (|^) we obtain a relation between the mean square 
spectra of the flows and the forcing in frequency space, 



I /I s 



uj 2 + (-Mm + g(k r )) 2 



(4) 



Assuming that |/| 2 is independent of k, u> (white noise), 
the integration of (||) over u) yields the relation between 
the mean square flow amplitude at an instant of time and 
the forcing, 



g{K)\ 



(5) 



The flow intensity (0) replaces the Bose distribution in 
the BEC case. Both functions tend to infinity when the 
amplification (stimulated emission) terms approach the 
damping (absorption) terms. As long as every mode is 
net-damped at a rate independent of the system size, the 
energy density stored in (0,0) modes must decrease pro- 
portional to the system size, since the total shear flow 
energy is distributed among an increasingly dense set of 
modes. However, analogous to the thermodynamic the- 
ory of the BEC, when the continuous flow spectrum is 
unable to hold the shear flow energy for non-zero min- 
imum net-damping rate and given random forcing, the 
nonlinear flow amplification term must adjust so that 
the remaining part of the flow energy is excited in the 
form of the most weakly damped modes, which are (0,0) 
modes. Hence, to demonstrate the possibility of conden- 
sation, it has to be shown that the flow amplitude in 
to ^ modes stays finite when the net-damping rate of 
the m = modes tends to zero, in the limit of infinite 
system size or, equivalently, in the approximation of a 
continuous poloidal mode spectrum. 

It is sufficiently general, to assume that g(k r ) has a 
maximum at k r = ko of order of the turbulence wavenum- 
bers and is parabolic near that maximum, g(k r ) — go — 
gi(k r — ko) 2 . The amplification terms will nearly can- 
cel the damping terms only for wavenumbers near fco, 
which justifies the approximation k r w fco which was 
made in the derivation of (|^). For the following anal- 
ysis we shift the fc r spectrum of the flows so that fco = 0. 
With ik r = d r , the operator in the denominator of (|^), 



(7 - .9o) - gidr + {ma {r ~ r mn )Y 



(6) 
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is the quantum mechanical Hamiltonian of the harmonic 
oscillator. Its eigenvalue for a mode with the "quantum 
numbers" (to, r mn , I), I S {0, 1,2,.. .}, is 

ui = 7 - 5o + 2^5i \aom\ (I + 1/2). 

The sum over I of the eigenmode contributions to (||) at 
fixed (m,r mn ) results in a logarithmic divergence, which 
stems from the infinitely broad random forcing spectrum 
and infinitely fast turbulence response. Hence, we cut 
off uji at an appropriate uj c depending on the turbulence. 
The sum is then approximated by an integral over I. The 
resulting amplitude associated with each pair (to, r mn ) is 



\v\ 2 (m,r mn ) 



\f\ 2 * 



In- 
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(7) 



with wo = 7 — go + ^fg~i\ctorn\ < uj c . The density of 
rational surfaces is |ao m l f° r given m. Approximating 
the sum over all \v\ 2 (m,r mn ) contributions to (jfy with 
to 7^ by an integral (which becomes exact for infinite 
system size), the total instantaneous energy density of 
the flow modes with m =/= 0, 



\a m\ ■ \v\ 2 (m,r mn )dm, 



(8) 



is obtained, where the integration interval is limited by 
the cutoff m c defined by ajo(m — to c ) = oj c . With the 
minimum net-damping rate = u>o(rn = 0) = 7 — go we 
obtain 



2|a |si 



- 1 + In 



(9) 



This expression converges to a finite value for Q — > 0. On 
the other hand, because the energy density of the m — 
modes, which tends to infinity for O — [the integral 
over k r of (||) does not exist for — iLUy ln + g = 0], has 
to be finite, we always have > 0. If the turbulence 
saturation requires a higher flow level than (^) at Q — ► 0, 
the description of the system by a continuum of poloidal 
mode numbers breaks down, the flow energy which can 
not be received by the to 7^ modes condenses into m = 
modes, and simultaneously r2 ^ 0. 

In a similar manner, it can be shown that in the limit of 
infinite system size the m = condensate is completely 
contained in the n = modes. Furthermore, there is 
no condensation of the radial wave numbers but the k r 
spectrum becomes arbitarily narrow around the point of 
weakest net-damping for large system size. 

Strictly speaking, condensation is unprovable by nu- 
merical studies, due to the restriction to finite system 
sizes. However, the validity of the individual parts of 
the model can be checked in the simulations. The local- 
ization of the m =/= shear flows on resonant surfaces is 
obvious in a plot of the flow spectra versus radius (Fig. |^) . 
The total flow amplitude associated with each (m,r mn ) 



quantum number in Eq. (Q) (Fig. ||) has much weaker 
slope than kg 2 oc mT 2 for kg,m — > 0. Therefore the in- 
tegral over all to 7^ shear flows (||) is expected to be 
finite in the limit of infinite system size (even if the esti- 
mate ^ for the individual amplitudes should be quanti- 
tatively wrong). Furthermore, the integral is reasonably 
well approximated by the corresponding sum in the finite 
system. Consequently, the infinite system will have ap- 
proximately the same ratio of to = flow amplitude to 
to 7^ flow amplitude. Finally, we note that the numer- 
ical studies agree with the above analytical prediction, 
that the flow condensate exhibits a strong peaking in k r 
for sufficiently large system size. 

Conclusions and consequences. - It has been shown 
numerically that in general the shear flows controlling 
the turbulence are not only (0, 0) modes but rather 
consist of a spectrum of poloidal mode numbers. The 
(to, n) ^ (0, 0) flows differ from drift waves or convec- 
tive cells by their large poloidal [10 times larger than the 
turbulence (Fig. S)] and parallel scale length, while their 
perpendicular scale length is similar to that of the tur- 
bulence. These shear flows are localized in the vicinity 
of resonant surfaces (Fig. ||) . In the limit of large system 
size, a non-zero (0,0)-mode amplitude develops only if 
the shear flows undergo a condensation into these modes, 
analogous to the Bose-Einstein condensation. Several 
features predicted by the analytic model have been re- 
produced by the numerical simulations. 

Due to the cancelling of damping and amplification 
terms, the (0, 0) flow condensate is practically undamped. 
This means in quantum mechanical language that the 
rate of absorption and incoherent re-emission, the "col- 
lision" rate, vanishes. Hence, far ranging interactions 
or ordering effects might be mediated via the shear flow 
condensate (but not by the uncondensed flows that suffer 
collisions and are pinned to resonant surfaces). As a con- 
sequence in the simple system used here the k r spectra 
become arbitrarily narrow. 

Since the flows depend on the distribution of ratio- 
nal surfaces and mode numbers, to accurately model 
the shear flow system in numerical studies, care has to 
be taken to not introduce spurious resonant surfaces or 
modes, e.g., by parallel extension of the flux tube [§,||. 
Remarkably, it can be shown that increasing the flux tube 
length does not lead to the correct limit of large system 
sizes, since, e.g., for an infinitely long flux tube, conden- 
sation into (0, 0) modes can not occur. 

Up to now, in flux tube based turbulence computa- 
tions the shear flows were implicitly assumed to be global 
modes. With domain widths too small for the large 
poloidal scales of the continuous part of the flow spec- 
trum, the flows appear to have zero poloidal and toroidal 
mode number. Such modes do not experience the resis- 
tive damping, which would reduce the flow amplitude in 
a full system. Hence, the simulations tend to overesti- 
mate the total flow amplitude, which may therefore ex- 
ert a strong stabilizing effect on the turbulence. To avoid 
an underestimate of the transport, flux tube simulations 
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have to be checked for influences of a finite poloidal scale 
length of the flows. 

The author would like to thank Dr. D. Biskamp for 
valuable discussions. 
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FIG. 2. Mean square shear flow amplitude as a function of 
kg for case (a) for Lg = 768Lo. Note the greatly different 
scale lengths of the turbulence (kg > 0.3) and the shear flows 
{k e < 0.3). 
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FIG. 3. Mean square shear flow amplitude as a function of 
kg and radius x = r — r for case (b) for Lg = 768L . Note 
the localization of kg ^ flows on resonant surfaces. The 
resonances lie on the intersections of families of hyperbolae in 
the x x kg plane, some examples of which are displayed. 
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FIG. 1. (0,0) shear flow energy density as a function of 
poloidal domain size Lg for case (a) (solid) without conden- 
sation and case (b) (dashed) exhibiting condensation; the thin 
line is proportional to 1/Lg. 
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FIG. 4. Mean square shear flow amplitude on a single 
resonant surface multiplied with ke oc m for case (b) for 
Lg = 768-Lo. The thin lines are proportional to kg 1 (steeper 
curve) and kg (flatter curve). 
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